function [f]=gen_f(nx,ny)
hx = 1/nx;
hy = 1/ny;
length = (nx+1)*(ny+1);
f = zeros(length,1);

for n=1:length
    [x,y] = ij(n,nx);
    f(n) = sin(pi*x*hx);
end
